Evolution of percolating force chains in compressed granular media 
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The evolution of effective force chains percolating through a compressed granular system is investi- 
gated. We performed experiments by compressing an ensemble of spherical particles in a cylindrical 
container monitoring the macroscopic constitutive behavior and the acoustic signals emitted by 
microscopic rearrangements of particles. As a novel approach, we applied the continuous damage 
model of fiber bundles to describe the evolution of the array of force chains during the loading 
process. The model provides a nonlinear constitutive behavior in good quantitative agreement with 
the experimental results. For a system of hard particles the model predicts a universal power law 
divergence of stress when approaching a critical deformation. The amplitude distribution of acoustic 
signals was found experimentally to follow a power law with exponent S = 1.15 ± 0.05 which is in a 
good agreement with the analytic solution of the model. 



> 

o 
o 
o 



X3 
O 

o 



X 



Recently, the behavior of granular materials has been 
extensively studied under various conditions due to their 
scientific and technological importance. Huge experimen- 
tal and theoretical efforts have been devoted to obtain a 
better understanding of the global behavior of granular 
media in terms of microscopic phenomena which occur at 
the level of discrete particles Subjecting a confined 

granular packing to an uniaxial compression a rather pe- 
culiar constitutive behavior can be observed: for small 
strains a strong deviation from the linear elastic response 
can be found implying that the system drastically hard- 
ens in this regime Linear elastic behavior can only 
be achieved asymptotically at larger deformations when 
the system gets highly compacted. When the external 
load is decreased again the system shows an irreversible 
increase in its effective stiffness, furthermore, under cyclic 
loading hysteretic behavior is obtained. 

Microscopically, inside a compressed granular packing, 
stresses are transferred by the contact of particles. Un- 
der gradual loading conditions the particles get slightly 
displaced changing their contacts and the local load sup- 
ported by them, which can be experimentally visualized 
using photoelastic materials for grains. These experi- 
ments revealed that in a compressed granular system the 
stresses are transmitted along the direction of the exter- 
nal load by force chains which can branch at the grains 
and form a complex network |Q| . Particles lying between 
lines of the force network do not support any load and 
can even be removed from the packing without changing 
its mechanical properties. Increasing the external load, 
more and more force lines appear and they all undergo 
erratic changes until the system reaches a saturated state 
when all the particles hold typically the same load and 
the system behaves as a bulk material. The creation 
and restructuring of percolating force chains implies rel- 
ative displacements of particles which can be followed 



experimentally by recording the acoustic waves emitted, 
however, up to now no such experiments have been per- 
formed systematically. Theoretically, this problem has 
been mainly studied by means contact dynamics simula- 
tions using spherical or cylindrical particles, and cellular 
automata Computer simulations also revealed the 

generation and evolution of force chains in compressed 
granular materials, however, the statistics of microscopic 
restructuring events, the emergence of the array of force 
chains and their relation to the macroscopic constitutive 
behavior remained un clarified. 




FIG. 1. Experimental set up and sketch of the array of 
force chains used in the model. 
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In the present letter the generation and evolution of 
percolating force chains is studied experimentally and 
theoretically in granular packings subjected to an uni- 
axial external load. We measured experimentally the 
macroscopic constitutive behavior and the acoustic sig- 
nals emitted by microscopic restructuring events com- 
pressing an ensemble of spherical glass bead confined in 
a cylinder. Based on the analogy of force lines percolat- 
ing through the system and fibers of a fiber composite 
we propose a novel theoretical approach, namely, an in- 
version of the Continuous Damage Model (CDM) of fiber 
bundles to describe the stress transmission through 
granular assemblies. The model naturally captures the 
emergence and gradual hardening of force chains and pro- 
vides analytic solutions for the constitutive behavior and 
acoustic activity as well. We also propose an efhcient 
Monte Carlo simulation technique. 
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FIG. 2. Constitutive behavior measured experimentally. A 
power law with exponent 2.6 was fitted to the data. The 
comparison with the analytic solution of the model is also 
presented. 



In our experiments a cylindrical container made of 
PMMA was filled with glass beads of 5 mm diameter and 
water. The cylinder has a thickness of 5 mm and a diam- 
eter of 140 mm. A punch test was carried out applying 
monotonically increasing displacements at the top level 
of the glass beads. Eight acoustic sensors were placed at 
the container wall to record the signals emitted durinj 
the compression of the beads, as can be seen in Fig 
The actual force and the displacement, measured at the 
traverse of the loading machine, and the acoustic signals 
were recorded simultaneously. An eight channel transient 
recorder was used as an analogue-digital converter to en- 
able the storing of the acoustic emission waveforms and 
a signal-based data. Experiments were performed un- 
der strain controlled conditions at a fixed strain rate, i.e. 
moving the traverse at a constant speed of 1 mm/minute. 

The nonlinear elastic response of the system can be ob- 
served in Fig. H where the measured force F is presented 
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as a function of strain e imposed. A power law of an 
exponent 2.6 was obtained as a best fit to the measured 
data in reasonable agreement with former experiments of 
Ref. §. 

The emergence and gradual hardening of force chains 
is responsible on the microscopic level for the strong 
non-linearity observed macroscopically. To obtain infor- 
mation about microscopic processes, the acoustic waves 
emitted due to sudden relative displacements of parti- 
cles were monitored. Typically several hundred signals 
were recorded during the experiment. The inset of Fig. 

shows the automatically extracted peak amplitudes of 
the burst signals versus time. The energy is defined as 
the integral of the acoustic emission signal amplitude fol- 
lowing the onset time. The energy values of the acoustic 
emissions are summed up in intervals of 30 seconds to elu- 
cidate the time dependent evaluation of acoustic emission 
activity. More details of acoustic emission data analysis 
and especially signal-based techniques can be found in 
The statistics of restructuring events is charac- 
terized by the distribution D{s) of the height of peaks s, 
which is presented in Fig. ||on a double logarithmic plot. 
It can be seen that D{s) shows a power law behavior 
over two orders of magnitude, the exponent of the fitted 
straight line is 5 = 1.15 ± 0.05. 
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FIG. 3. The statistics of acoustic signals. A power law of 
an exponent 5 — 1.15 ±0.05 was fitted to the size distribution 
of the signals of the inset. 

Based on the analogy of force chains in granular pack- 
ings and fibers of fibrous materials or fiber reinforced 
composites, we propose a novel theoretical approach, 
namely, an inverted fiber bundle model, to describe the 
gradual emergence and hardening of force chains occur- 
ring under uniaxial compression. Fiber bundles are com- 
posed of parallel fibers of identical elastic properties but 
stochastically distributed breaking thresholds. A fiber 
fails during the loading process when the local load on it 
exceeds its breaking threshold. Fiber failures are followed 
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by a redistribution of load on the remaining intact fibers 
according to the range of interaction in the system. The 
so-called Continuous Damage Model introduced recently 
1^,^, is particularly suited to model granular materials 
since it captures gradual stiffness changes of elements of 
the model. In our model of compressed packings, force 
lines formed by particles are represented by an array of 
lines organized in a square lattice as illustrated in the 
inset of Fig. |^. A randomly distributed rearrangement 
thresholds d is assigned to each line of the array from a 
cumulative probability distribution Pq{-^), where de- 
notes the characteristic strength of force lines. During 
the compression process, when the local load on a line 
exceeds its threshold value d the line undergoes a sudden 
restructuring as a result of which it becomes stiffer and 
straighter. The lines' stiffness increases in a multiplica- 
tive manner, i.e. the stiffness is multiplied by a factor 
a > 1 at each restructuring so that the constitutive equa- 
tion of a single line after suffering k restructurings reads 
as fj = Eott^e". Here Eq denotes the stiffness character- 
izing single particle contacts, and the exponent a takes 
into account possible non-linearities of a single contact 
like for Hertz law a = 1.5. After each rearrangement 
the force chain gets a new threshold value (annealed dis- 
order) from a distribution of the same functional form, 
but the characteristic strength dc of the distribution is 
increased in a multiplicative way so that after k rear- 
rangement events the disorder distribution takes the form 
= ^0 ( ) ■ The maximum value of possible re- 
structurings kmax is proportional to the number of con- 
tacts and therefore to j- in the nomenclature of the inset 
of Fig. The ratio r = a/g is a very important pa- 
rameter of the model, it decides whether the force chain 
becomes more fragile (r > 1) or more ductile (r < 1) 
as a result of restructuring. Essentially, t is material 
dependent and determined by the stability of contacts. 
Calculations showed that the most interesting behavior 
of the model can be found when r = 1, hence, in the 
present letter our investigation is restricted to this case. 
The full complexity of the model will be explored in a 
forthcoming publication [T^ ]. 



Experiments and discrete element simulations 




have revealed that the number of effective force chains 
increases during the compression process until it reaches 
a saturation value. To capture this effect in our model, 
for the number of elements we prescribe the form N{e) = 
NoG{e), where No denotes the saturation number of 
chains, and the profile G{e) has the property G{e) — > 1 
with increasing e. Hence, the number of force lines 
dN emerging due to an infinitesimal deformation in- 
crement from e to e -I- de is dN = Nog{s)de, where 
g{e) = dG{e)/de. Following the derivation of the consti- 
tutive behavior of the continuous damage model of fiber 
bundles [^|j9|, the macroscopic constitutive equation of 
the compressed granular system can be cast into the form 



k — 1 

a{e)^E/Y. f a^{e-e*rg{e*)P^^{e-e*)x (1) 
[l-Po{e~e*)]de* 

+ r a'^— (e - e*)"(7(e*)Po*'""^(e - e*)de* . 
Jo 

Eq. (|l|) takes also into account that the local strain of 
force lines e — s* is different from the externally im- 
posed strain value e since it also depends on the initial 
strain e* . The integral is performed over the whole load- 
ing history to take into account all the lines generated. 
For explicit calculations we imposed an exponential form 
N{£) — No{l — e~^) for the number of chains. The best 
fit obtained to the experimental data is presented in Fig. 
H where the force F = NqCt is plotted against deforma- 
tion e. It can be seen that a good fit was achieved with 
physically reasonable parameter values. The maximum 
possible number of percolating force chains No that can 
emerge in the system was estimated as the ratio of the 
total area of the container Ao to the cross section of a 
single particle Ap, i. e. No — ^ = 784. The value of the 

other parameters are Eo = 4600;^, do = 3A^, /? = 0.01 
and a = q = 1.01. The value of a falls close to one in- 
dicating that a single restructuring gives rise only to a 
slight increase of stiffness of a force chain. Model calcu- 
lations revealed that the zero derivative at the starting 
part of the constitutive curve is due to the gradual cre- 
ation of load bearing force chains. The small value of 
(3 implies that the generation of new force chains stops 
at a relatively small strain value, and hence, the later 
rapid increase of F as function of e is mainly caused by 
the hardening of the existing force lines. Further infor- 
mation can be gained about the constitutive behavior 
for larger strains by simplifying Eq. (|^) assuming a fixed 
number No of force chains from the beginning of the pro- 
cess. Under this assumption Eq. (|l]) can be reformulated 
as 

a{e)=Eoe[l-Po{e)] £ a'PoisY (2) 

i=0 

It can be seen from Eq. (|^) that if the maximum number 
kmax of possible restructuring events goes to infinity the 
stress cr has finite values only for aPo{e) < 1. In this 
case the summation can be performed in the first term, 
while the second term tends to zero, and the constitutive 
equation takes the form 

.(s)=i..s[l-Po(e)]^-^. (3) 

It follows that the stress a diverges when e approaches 
a critical value Sc, where Sc satisfies the equation 
Poisc) — l/a. Expanding Po{e) into a Taylor series at 
£c as Po(e) = Po{£c) +J'o(£c)(ec — e) + . . ., and substitut- 
ing it into Eq. (^) the behavior of a in the vicinity of Sc 
reads as 
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a(e) 



1 



(4) 



It means that the stress a shows a power law divergence 
when e approaches the critical value Sc- The value of the 
exponent is universal; it does not depend on the form of 
disorder distribution Pq, while the value of Ec depends 
on it. It is interesting to note that in Rcf. |^ the same 
power law divergence was found in large scale molecular 
dynamic simulations of a hard sphere system. 

When the number of force lines is fixed it is possible to 
obtain analytic results also for the statistics of restruc- 
turing events. Restructuring occurs during the compres- 
sion process when the local load on a force line exceeds 
its threshold value. Since loading is performed under 
strain controlled conditions, there is no load redistribu- 
tion among existing force lines, i.e. restructuring of a 
force line does not affect other elements of the system. If 
the new threshold value, assigned to the force line after 
rearrangement, is smaller than the local load, the force 
line undergoes successive restructurings until it gets sta- 
bilized. The number of steps to reach the stable state 
defines the size s of the restructuring event, which is the 
analog of the acoustic signals measured experimentally. 
The number of restructuring events Nk^s of size s starting 
in force chains which have already suffered k restructur- 
ings can be deduced as 



No 



= Po{e)Pr'-'ie) [l-Po{e)], 



for s + k < kmax — 1, and 



(5) 



(6) 



for s + k ~ kmax (see also Ref. |^] ) . The number of events 
D(s) of size s can be determined by integrating over the 
entire loading history and summing over all possible k 
values 
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No 



de. (7) 



Finally, substituting Eqs. (^,^ into Eq. and perform- 
ing the calculations yields 
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where 



1 < s < fc„ 



(8) 



i.e. the distribution of microscopic restructuring events 
exhibits an universal power law behavior with an expo- 
nent 1, which is completely independent on the disor- 
der distribution. Numerical simulations revealed that the 
universal power law behavior also holds when the gradual 
creation of force chains is taken into account, i.e. when 
the system is described by the full Eq. (^. The statistics 
of restructuring events obtained by Monte Carlo simula- 
tions is presented in Fig. ^ where simulations were per- 
formed under the assumption N{e) = iVo(l — ). The 



inset shows local events of different sizes that occurred 
during the loading process, and their distribution is pre- 
sented in the main figure. The power law behavior of the 
analytic prediction of eq. ^ is verified. It is important to 
emphasize that the theoretical results on event statistics 
(Fig. ^ are in a very good quantitative agreement with 
the experimental findings (Fig. ^). 
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FIG. 4. Inset: local events of size s from a Monte Carlo sim- 
ulation of strain controlled loading. The statistics of events 
D{s) is presented in a double logarithmic plot. 

Note, that the functional form and the value of the 
exponent of D{s) in the analytic calculations is mainly 
the consequence of the locality of restructurings due to 
the absence of load redistribution. The excellent agree- 
ment observed indicates that this is likely the microscopic 
mechanism responsible for the power law statistics of 
acoustic signals observed experimentally. 
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